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A global search for possible LiF cluster structures is performed, up to (LiF)g. The method is based 
on simulated annealing, where all the energies are evaluated on the ab initio level. In addition, the 
threshold algorithm is employed to determine the energy barriers for the transitions among these 
structures, for the cluster (LiF)4, again on the ab initio level; and the corresponding tree graph is 
obtained. 



I. INTRODUCTION 



-^ . The question of the possible structures of a chemical system, both in the case of molecules and of solids, is of huge 

importance [lH8( ■ Identifying these structures requires the determination of the local minima on the energy landscape 
of the system, or more generally the locally ergodic regions (basins). At low temperatures, individual local minima, 

r~ — i- that are surrounded by sufficiently high energy barriers, are kinetically stable. There are numerous studies of clusters, 
where global searches for their minima have been performed. This is a non-trivial task, as the number of minima is 
^ . expected to increase exponentially with system size, see e.g. [9l.[l0|. 

r—i ' However, finding the local minima is only one part of identifying promising structures. The second task is the 

^— i , analysis of the barrier structure of the landscape, in particular the determination of the energetic barriers separating 
the local minima, which control the stability of the structure. To describe the energy landscape and barriers of a 
system, tree graphs [8L H1H13J can be employed. Examples for such tree graphs for solids can be found in [12, H]) and 
Cd ' for clusters in @, respectively. 

In the present work, we have chosen LiF clusters as the system to be studied, since alkali halide clusters have served 
as important model systems for complex energy landscapes. Such clusters have previously been investigated with 
various methods. In early calculations, atoms were placed in random positions and then moved (to be more precise, 
a quench was performed, i.e. a move is suggested and accepted, if it lowers the energy) [15| . A different scheme 
applied for NaCl clusters [16| was to build up clusters by adding one atom at a time: an atom is added in a random 
location, then the energy of the new structure is computed, and according to an acceptance probability, this site is 
accepted or not. When the requested number of atoms has been added and the intended cluster size is obtained, then 
subsequently a gradient minimization of the cluster is performed. With a sufficiently high number of trials, this is 
supposed to give a good overview of the energy landscape. A further technique used was simulated annealing with 

\^\ | molecular dynamics or simulated annealing where the system is represented by a Gaussian density distribution for 
each particle [l7[- Basin hopping was applied [lj| to charged clusters (NaCl) ra Cl~, and to the related system Al ra N ra 

C*~) ' [lj| . In addition, a genetic algorithm has been applied to NaCl clusters [2fJ, [21| . 

r — . The aforementioned studies employed model potentials. A comprehensive ab initio study on the level of Hartree- 

f^i ■ Fock, second order perturbation theory and the coupled pair functional has been performed for selected clusters 
[22l | . Moreover, ab initio molecular dynamics simulation was applied to small NaCl clusters (stoichiometric and non- 
stoichiometric) [23|. Ab initio calculations on known structures have nowadays become routine work, see e.g. [24JI25J. 
Among the studies mentioned so far, most apply to alkali halides in general, or to NaCl clusters as a special case. LiF 
clusters are addressed in references Y2M [26l. ut\ . For a review on the molecular structure of metal halides, see [28l[29|. 
and on the modeling of clusters, see [301]. 

The search for local minima is commonly split in two parts [l2j: a global search for structure candidates, and a 
subsequent local optimization with accurate ab initio energy calculations. Usually, the global search is performed 
with empirical potentials, because the global search is the time-consuming part, and empirical potentials require little 
CPU time. In contrast, in the present work, both steps, i.e. global search and local optimization are performed on 
the ab initio level. The procedure is still split in two steps: first, the global optimization is performed with ab initio 
calculations on a lower level of accuracy which makes them much faster, however. Nevertheless, this is a non-trivial 
task as the CPU times are still much higher than with empirical potentials. The local optimization is then performed 
with the usual high accuracy, as this step requires less CPU time than the global optimization. Employing ab initio 
methods in both steps of the search has several advantages: the most important is that ab initio methods work for 
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general systems, whereas empirical potentials work reasonably well for ionic systems, but less so for covalent systems; 
and some knowledge of the expected bond type is required in advance to choose the potential. 

Recently, we demonstrated, in the context of the structure prediction of solids, that a full ab initio treatment is 
feasible in both stages, i.e. the global search and the subsequent local optimization can both be performed on the 
ab initio level. LiF was chosen as an example for an ionic solid [3l| , and BN as an example for a covalent solid 32] . 
In the case of lithium fluoride, it was found that the relevant minima were the same when full ab initio structure 
prediction was performed, compared to the earlier studies with model potentials [33| (for a brief summary see also 
[34j ) . Concerning clusters, the earliest work which combined simulated annealing with ab initio energies employed the 
Car-Parrinello method [35|], and was performed for e.g. silicon [36| and selenium [37| clusters. In order to include the 
nuclear motion explicitly, ab initio path integral molecular dynamics |38l44C| was suggested as an extension. In these 
methods, the atom moves are determined by molecular dynamics. Simulated annealing in combination with stochastic 
moves and ab initio energies had been applied earlier to various clusters, e.g. the clusters LisH [41] . Mg + (H2 0)„ [42[ 
or lithium clusters [43j . 

The goal of the present article is the study of the energy landscape of LiF clusters, completely on the ab initio level. 
The first step is to study the feasibility of a global optimization with simulated annealing and employing ab initio 
energies for these systems, and to compare with earlier studies as far as available. The second and more challenging 
step is to determine the energy barriers between the structures, again on the ab initio level. The latter task is fairly 
time-consuming, and had not been attempted for the solids which had been studied previously with ab initio simulated 
annealing. The intention of the present article is to test and demonstrate the feasibility of the latter step. LiF has 
been chosen, since the binding situation is fairly ionic, and thus convergence for unknown structures is relatively 
easy to achieve. Thus, comparatively simple ab initio methods such as Hartree-Fock or density functional theory are 
expected to work reasonably well. The outcome of the present study, and also the information about the demand 
on computational resources, is important in order to estimate the effort for more difficult systems with a partially 
covalent character, which will be future targets. 

In the following paragraph, the technical details will be described. Then, the results are presented and discussed; 
and finally the article will be summarized. 

II. METHOD 

As mentioned above, the general optimization procedure consists of several steps: first, a simulated annealing run 
with a subsequent stochastic quench is performed, to identify possible candidate structures. This is followed by a 
local optimization based on analytical gradients. This is repeated many times, in order to identify as large as possible 
a set of structure candidates, and to obtain some statistics about the structures found. 

The details for the LiF calculations are described in the following paragraphs. Between 1 and 8 formula units of 
LiF, i.e. from 1 lithium and 1 fluorine atom up to 8 lithium and 8 fluorine atoms were placed at random positions in 
a large box. The box was initially cubic with a size corresponding to three times the total volume of the atoms/ions 
as estimated from the atomic/ionic radii. Note that no periodicity is used, i.e. the box is not repeated in space. The 
box is mainly used to place the atoms initially, and, in general, moves are defined which allow to change the box size. 
In the case of molecules, one might consider not allowing such moves and keep the size of the box fixed, which might, 
however, lead to the search requiring more steps to reach a local minimum. 

The length of the simulated annealing run was in the range between 5000, for the smaller clusters, and 75000 steps, 
for the largest clusters, respectively. The initial temperature was chosen in the range of 1-10 eV, (1 eV corresponds 
to 11604 Kelvin), and was reduced by up to ~ 50% during the simulation, in the longer runs. This temperature is 
used when performing the Metropolis Monte-Carlo simulation: the energy difference per atom AE is divided by this 
temperature when evaluating the expression exp(— AE/ksT). 

The simulated annealing was followed by a quench with 10000 steps, i.e. a simulated annealing run with a temper- 
ature of eV, which means that only downhill moves are allowed during the quench. 

The moves were chosen as: moving individual atoms (70%), exchanging atoms (10%), changing the size of the 
simulation box (20%). Note that exchanging atoms is a useful move especially if the initial random structure has 
neighboring pairs of the type Li-Li or F-F. No symmetry was prescribed during the simulated annealing and quench 
runs, i.e. the point group was always G\. 

A minimum distance between two atoms (given by the sum of the radii of the atoms, multiplied by 0.7) was 
prescribed in order to avoid unrealistic geometries which may lead to numerical instabilities. The radii used were 
based on tabulated values for atomic and ionic radii, as a function of charge, and the Mulliken charge computed for 
the previous configuration. In those moves which change the size of the simulation box, the probability of reducing 
the lattice constant was enlarged to 70%, to speed up the reduction of the cell size. 



With this choice of the parameters, one simulated annealing and quench run lasted from minutes for the smaller 
clusters up to 1 day for the largest cluster considered, i.e. (LiF)g. 

The ab initio calculations were performed with the CRYSTAL06 code [44|], which is based on local Gaussian type 
orbitals. The basis sets employed are displayed in table |U Two slightly different basis sets are used during the 
global search and the local optimization. During the global search, slightly tighter sp functions were chosen for the 
two outermost fluorine exponents, in order to enhance the numerical stability and the speed of the calculations. 
In order to significantly reduce the CPU time, the following simplifications were used during the global search: 



the thresholds for the integral selection |44j were reduced from the default values (10 6 , 10 6 , 10 6 , 10 6 , 10 12 ) to 
10~ 4 , 10 -4 , 10~ 4 , 10~ 4 , 10 -8 , and the convergence threshold of the self-consistent field cycles was reduced from 10~ 7 
(default) to 10~ 3 . 

The global search was performed on the level of Hartree-Fock theory, i.e. all the energies of all the structures 
appearing during the simulated annealing and quench were computed from first principles. As has been demonstrated 
earlier [323, Hartree-Fock theory has the advantage that convergence for random structures is facilitated due to the 
large HOMO-LUMO gaps (gaps between highest occupied and lowest unoccupied molecular orbital) on this level of 
theory, compared to the much smaller LDA (local density approximation) gaps. 

The local optimization is essentially a routine task and employed analytical gradients as implemented in the CRYS- 
TAL06 release |4a - l50l ]. This local optimization was performed on the level of Hartree-Fock, and subsequently, the 
optimal Hartree-Fock structure was again optimized on the LDA level. In very few cases, a structure found on the 
Hartree-Fock level transformed to a different structure on the LDA level. These structures were dropped, and only 
those structures were included which remained the same after the LDA optimization. This way, a consistent set of 
local minima was obtained for clusters of the size (LiF)i to (LiF)g. The symmetry of these clusters was analyzed with 
the program SYMMOL [H|. 

Going beyond the search for local minima, an exploration of the barrier structure using threshold runs [12j | was 
performed for the system (LiF)4. A threshold run starts from a local minimum and explores the part of the energy 
landscape that can be reached from this configuration without crossing a given energy lid. A new configuration 
is generated (according to a certain move class), and this configuration is accepted if its energy is below a certain 
threshold. Then, one or several quenches are performed starting from some of the high energy structures encountered 
during the threshold run. When low thresholds are applied, these quenches tend to return to the starting minimum. 
However, they can reach a new structure, if the threshold energy is high enough to overcome the barrier between the 
initial structure and the new one. This way, an upper bound of the energy barrier between two structures can be 
determined. 

Here, five threshold runs with a length of 100000 steps each were performed for each energy threshold and each 
starting minimum. The runs were interrupted after 20000, 40000, 60000, 80000 steps, in order to perform three 
quenches with 20000 steps each. This results in a total of 400000 steps per threshold run (100000+5*3*20000). Note 
that a quench is again a stochastic procedure, and thus three quenches can, in principle, end up in three different 
local minima. This, however, turned out to be only rarely the case (in about 1% of the cases; in all the other cases, 
the same minimum structure was obtained from the three quenches) . The thresholds were chosen in an iterative way 
so that the energy barriers could be determined to ~ 0.01 eV/atom. A more systematic approach employing a set 
of threshold values as in [14| would also provide more statistical information such as the probability to find a certain 
structure as a function of the threshold energy, but would require significantly more CPU time. 

As the threshold run intends to give a realistic simulation of the transition from one cluster structure to another 
one, the exchange of atoms was not permitted, and instead, only atom moves (80%) and changes of the box size (20%) 
were allowed. 

With the present choice of parameters (employing the set of less accurate parameters, as used in the global search), 
the CPU time for one energy calculation of a (LiF)4 cluster is about 1 second on a single CPU (Intel Xeon 5150, 2.66 
GHz). More precisely: in the simulated annealing run, the initial geometry is like a gas, the atoms are far apart, 
and many integrals are discarded (due to the selection criteria of the integrals). Thus, one energy calculation at this 
geometry takes only a fraction of a second (less than 0.2 seconds). In the final geometry at the end of the simulated 
annealing run, the atoms are close, and many more integrals have to be evaluated: as the atoms are closer, the overlap 
integrals are larger, and fewer integrals can be neglected. One energy calculation takes now up to 2 seconds. During 
the threshold run, the CPU time for one energy calculation is of the order of 1-2 seconds. 

One threshold run with 400000 energy calculations thus lasts of the order of 400000 seconds (5 days), and thus, 
with 5 threshold runs per lid, about 4 weeks of CPU time are required. This can obviously be trivially parallalized by 
employing 5 CPUs, one CPU for each threshold run. An initial threshold may be estimated as e.g. 0.5 eV/atom, and 
subsequently this value can be bisected, in order to determine the energy at which a transition becomes feasible. The 
bisection is repeated, until the desired accuracy of the barrier is achieved. The information is insofar redundant, as a 
transition between two cluster types can be observed in both directions: e.g. when the barrier from the cube-shaped 
cluster to the ring-like cluster is known, then this can be used to have an initial guess for the threshold energy which 



is to be used for the ring-like cluster. 

Similarly, the guess of the initial threshold can be better adjusted when some results from other clusters become 
available (i.e. when the order of magnitude of the barrier is known). In total, for each of the 5 different (LiF)4 
clusters, about 5 threshold values had to be selected, in order to determine the tree graph. This results in a total of 
25*4 weeks of CPU time on a single CPU for the whole tree graph. 

III. RESULTS AND DISCUSSION 

A. Structures found 

The most relevant structures found for the (LiF) ra (n = 1,...,8) clusters are displayed in figures 1-6, visualized with 
XCrysDen [52J . The statistics are compiled in table [TTJ Concerning comparisons with other results for alkali halide 
clusters, it should be noted that most of the available literature data refers to NaCl clusters, and there is only few 
data available explicitly for LiF clusters [24l l26l l27l|. However, due to homology, one would expect similar minimum 
structures to be present in LiF and e.g. in NaCl clusters. 

For n = 1, 2 and 3 formula units, only one structure was found for each system (figure [TJ that constituted a stable 
minimum for both HF (Hartree-Fock) and LDA calculations. Additional modifications that were stable only on the 
HF level appeared one time each for n = 2 and 3. These structures were thus dismissed. A double-chain like structure 
was never observed. This confirms the statement in [2J] that such a double chain is stable for Na3Cl3 [22|, but not 
for Li 3 F 3 . 

For four formula units, four structures were found in 75 simulated annealing and quench runs, and a fifth during 
the threshold runs. These structures agree with those given for (NaCl)4 clusters in [22] (apart from the structure in 
figure lc of [22] which was not observed). 

For n = 5, seven energetically low-lying structures were found (figure [3]). Three of them correspond to the set of 
minima of (NaCl) 5 as published in references [16l.l22j]; the lowest structure in energy agrees with the global minimum 
of (NaCl)s: a cube with one additional LiF attached to one edge. We note that one particularly low-lying structure 
(5c) was found which had not been observed in the earlier work on NaCl clusters (l6j |. 

Figure [4] presents the eight structures with lowest energy obtained for (LiF) 6- Structure 6a with two six-membered 
rings atop of each other is the energetically most favorable one. This structure is a fragment from the 5-5 [33| (or 
'hexagonal MgO') structure which is a hypothetical low-lying modification of LiF [31| . This result is very interesting for 
experimental chemistry, since it suggests that the 5-5-type modification of LiF might be accessible via the deposition 
and growth of selected (LiF) n -clusters on an appropriate substrate (see also [53j). The structure with the second 
lowest energy is a cuboid cut from the rock salt structure which is the modification observed for bulk LiF at standard 
conditions. Structure 6c consisting of a cube, with a square attached, is ranked number three in energy. Note that 
in [161 |22j, a similar structure is found where the 3x2 rectangle in 6c is opened up. This structure was also found 
in the present work and is a candidate structure on the HF level. When optimizing it on the LDA level, it became 
the structure 6c, however. Thus, this additional structure was not included in the statistics. Besides those structures 
in table IIII1 additional minima with higher energies were observed. Furthermore, one would expect more structures 
resembling 6d to exist, i.e. consisting of a cube with two LiF units attached to different edges of the cube. Finally, a 
twelve-membered ring structure is now energetically less favorable than the more compact structures, which is to be 
expected. 

For n = 7, the lowest energy structure found agrees with the one from the literature for (NaCl)7 ]16,[2J] and (LiF)7 
[2J, [26| . In figure [5J the eight structures found with the lowest energy are shown. Note that some of the structures 
exhibit chirality. In some cases (7c and 7d, 7f and 7g), both enantiomers were found in the simulations, though not 
in the case of the structure 7b. The latter minimum was found only once, and we expect that increasing the number 
of runs should also produce the second enantiomer. Chirality can thus be used as a test of how thorough the energy 
landscape has been searched, since one would expect to find both enantiomers with the same probability. 

The largest clusters investigated contained eight LiF formula units. The lowest in energy are displayed in figure [U 
The energetically most favorable structure (8a) has S4 symmetry, and corresponds to the third lowest energy isomer 
of (NaCl)8 [22[ . The cuboid is next in energy and was lowest in energy for (NaCl)g [22|]. However, two energetically 
low-lying minima listed in earlier studies were not observed in the present work: an eight-membered double ring with 
D^ symmetry (for a corresponding NaCl cluster, see [22j, and for LiF, see [24], [26|, l27j). and a structure made of 
three cubes, with C s symmetry (for a corresponding NaCl cluster, see [22| . and for LiF, see [2g). To check their 
importance, the energies of these minima were computed, yielding -856.2474 Eh (HF) and -853.8478 Eh (LDA) for the 
structure with D\d symmetry and -856.2275 Eh (HF) and -853.8433 Eh (LDA) for the structure with C s symmetry, 
respectively. This would put the latter two structures in the second and fourth place (HF), or third and fourth 



place (LDA), respectively. It is expected that these minima could be found if more simulated annealing runs were 
performed. 

In figure [7] and table IIII| the energies (per formula unit) of the most favorable structures are compared. They 
decrease monotonously with cluster size and thus it is always energetically more favorable to form a bigger cluster 
from two smaller ones. However, if one considers the difference in chemical potentials (approximated by energy 
differences fx n = E n — E n -i, A^„ = fj, n+ i - fj, n = (E n+1 - E n ) - (E„ - E n -x) = E n +i - 2E n + E n -i), one notes 
that for certain cluster sizes it is slightly favorable to combine a (LiF) n +i and a (LiF)„_i cluster to form two (LiF) n 
clusters. This is also displayed in figure [7J A.fj, n is positive for n = 4, 6 and negative for n = 5, 7 which indicates that 
the clusters with n = 4, 6 have a higher stability. The reason appears to be that structures which resemble a fragment 
from the bulk are more favorable, i.e. cuboids or distorted cuboids, which cannot be constructed with 5 or 7 formula 
units (for a discussion of this issue, see also [a, [28|). This is reminiscent of magic cluster numbers in metallic and 
intermetallic clusters, but due to the ionic nature of the Li-F bonds, the effect is less pronounced. 

We note that with increasing size, rings become less favorable compared to the other configurations which are more 
bulk-like. This is to be expected, as their energy will approximate the energy of a one-dimensional chain. 

Concerning the statistics, for the larger clusters (i.e. n > 5), the percentage of successful runs becomes smaller; i.e. 
those runs which end up in one of the lowest energy minima in table [II] Quite generally, the structures of the global 
minima found usually agree with the ones known from the literature for LiF, or for NaCl. The only exception is 
(LiF)s, where the minimum found in the present work is more favorable than the one previously suggested for (LiF)g 
Ellillil. 

B. Energy barriers 



The barrier landscape of the system (LiF)4 was investigated with the threshold algorithm [12(, on the Hartree-Fock 
level. Note that all these calculations were performed with the parameters used during the global search, as these 
calculations are very time-consuming, and the weaker parameters used during the global search lead to a significant 
speed-up. No local optimizations were performed, for reasons of consist ency : during the threshold run, the energies 
were computed with less accurate ab initio parameters (as explained in |3l|), and thus also the energies of the local 
minima were computed on this level. The energies displayed in figure [8] are thus slightly different from the ones in 
table HIl From the results of the threshold runs, the tree graph for the ab initio energy landscape of (LiF)4 was 
constructed (see figure [SJ. 

We note that the cube (4a) and the ring-like structure (4b) have a relatively high stability and the barriers are 
larger than 0.1 eV/atom (i.e. 0.8 eV in total). A crude estimate of the lifetime of such a cluster can be performed by 
applying Arrhenius law k = Aexp(—EB/kBT): if we assume a prefactor of A — 10 13 i (in the range of the vibrational 
frequencies, see j54[), then at a temperature of 300 K, a barrier Eb of 0.8 eV would result in a rate constant k of 
the order of 0.4 -. Thus, one would estimate that both structures might be observable in the experiment. There is 



indeed some experimental evidence for a ring-like (LiF) 4 cluster 54|. The other three structures found have virtually 
no energy barrier, and easily transform to the ring or cube structure. This is to be expected, as one can see that the 
rectangular structure (4c) can easily be deformed to become ring- like (4b). In addition, the two very high lying minima 
(4d and 4e) were so rarely observed during the global search and the threshold runs, that the basins corresponding to 
these structures have either a low energy barrier or a very small volume in configuration space. This is consistent with 
the following fact: when threshold runs were started in structure 4d or 4e, with an energy lid only slightly above the 
energy of these minima, then the walkers only rarely returned to the starting minimum, but instead usually changed 
to the ring, cube or rectangular structure. 

IV. CONCLUSION 

Employing simulated annealing as a global optimization method, and using ab-initio energies during both the global 
and the local optimization, the structures of lithium fluoride clusters (LiF)„ (n — 1,...,8) have been predicted. In 
particular, the lowest energy structure agreed for all cluster sizes with, or was more favorable than, the known (LiF)„ 
minimum structures. The majority of the low-lying minima previously known from calculations on various alkali 
halide clusters such as (NaCl)„ was found. 

As a new methodological development, the threshold algorithm was applied to a molecular system, (LiF)4, where 
all the energy calculations were performed on the ab initio level. This way, the energy barriers separating the minima 
of the energy landscape of this cluster were computed, and a tree graph representation of the energy landscape of 
LiF4 was constructed. Based on this tree graph, two of the cluster modifications may be accessible in the experiment. 



Estimating the stability of the various clusters raises the general issue of the accuracy of the calculated barriers 
separating the cluster modifications. Clearly, in the ideal case, one would employ highly accurate quantum chemical 
methods, such as coupled-cluster calculations and possibly multi-reference methods. However, the determination 
of the barriers requires very large computational resources, and is thus not feasible on this higher level of theory. 
Thus, some compromise between accuracy and computational expense must be established. In the present study, 
the barriers have been computed on the Hartree-Fock level, where some additional calibration work was necessary, 
because standard Hartree-Fock calculations would be too slow. This type of calculation is orders of magnitude more 
expensive than employing an empirical potential, yet still not as accurate as coupled-cluster calculations (see, e.g. 
[551]). Nevertheless, the present work also demonstrates that the global exploration of the barrier structure of the 
energy landscape of small systems is feasible on the ab initio level, though not with the desirable highest accuracy. 
In particular for systems where no trustworthy empirical potentials are available that can well describe the barrier 
structure, such as e.g. compounds with a more covalent character, the present approach may well prove to be fruitful 
and significantly improve the quality of the energy landscape explorations, in spite of its current limitations. 
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TABLE I: Basis sets used for the global search (I) and the local optimization (II). 



basis set I 




basis set II 


exponent 


contraction 


exponent 


contraction 






Li 










s 






840.0 


0.00264 




840.0 


0.00264 


217.5 


0.00850 




217.5 


0.00850 


72.3 


0.0335 




72.3 


0.0335 


19.66 


0.1824 




19.66 


0.1824 


5.044 


0.6379 




5.044 


0.6379 


1.5 


1.0 


sp 


1.5 


1.0 


0.525 


1.0 1.0 




0.525 


1.0 1.0 






F 










s 






13770.0 


0.000877 




13770.0 


0.000877 


1590.0 


0.00915 




1590.0 


0.00915 


326.5 


0.0486 




326.5 


0.0486 


91.66 


0.1691 




91.66 


0.1691 


30.46 


0.3708 




30.46 


0.3708 


11.50 


0.4165 




11.50 


0.4165 


4.76 


0.1306 


sp 


4.76 


0.1306 


19.0 


-0.1094 0.1244 


19.0 


-0.1094 0.1244 


4.53 


-0.1289 0.5323 


4.53 


-0.1289 0.5323 


1.37 


1.0 1.0 


sp 


1.37 


1.0 1.0 


0.45 


1.0 1.0 


sp 


0.437 


1.0 1.0 


0.20 


1.0 1.0 




0.147 


1.0 1.0 



TABLE II: The total energies of the various clusters, in hartree units, and statistics of the outcome of the simulated annealing 
runs. A run is counted as successful, if one of the low lying minima was found. 



number of configu- 
formula ration 
units 



HF 



LDA 



-428.0742 
-428.0748 
-428.0616 
-428.0396 
-428.0336 



■426.8772 
-426.8608 
-426.8554 
■426.8384 
■426.8308 



success 



ful 



. / total ] 



-856.2524 -853.8572 



-856.2383 
-856.2216 
-856.2216 
-856.2199 



-853.8552 
■853.8226 
■853.8226 
■853.8181 



-856.2182 -853.8159 
-856.2142 -853.8152 



success 



ful 



. / total runs: 



number of 
times found 



1 


a 


-106.9444 -106.6406 


1/1 


2 


a 


-213.9925 -213.3883 


30/30 


3 


a 


-321.0393 -320.1310 


60/60 



41 
23 
10 
1 
° 
75 / 75 



5 


a 


-535.1039 


-533.6054 


15 




b 


-535.0961 


-533.5891 


1 




c 


-535.0979 


-533.5888 


15 




d 


-535.1058 


-533.5862 


7 




e 


-535.0815 


-533.5766 


2 




f 


-535.0859 


-533.5765 


1 




g 


-535.0725 


-533.5685 


1 




successful runs / total runs: 


42 / 50 


6 


a 


-642.1672 


-640.3706 


46 




b 


-642.1567 


-640.3670 


2 




c 


-642.1353 


-640.3364 


6 




d 


-642.1317 


-640.3316 


4 




e 


-642.1360 


-640.3238 


18 




f 


-642.1294 


-640.3221 


10 




g 


-642.1250 


-640.3172 


4 




h 


-642.1347 


-640.3101 


3 




successful runs / total runs: 


93/120 


7 


a 


-749.1967 


-747.1024 


6 




b b 


-749.1856 


-747.0938 


1 




c 


-749.1929 


-747.0932 


12 




d 


-749.1929 


-747.0932 


11 




e 


-749.1888 


-747.0920 


4 




f 


-749.1755 


-747.0753 


2 




g 


-749.1755 


-747.0753 


1 




h 


-749.1714 


-747.0700 


4 




successful runs / total runs: 


41 / 80 



16 
3 
3 

1 
1 
2 
1 
27 / 60 



a This structure was found only in the threshold run. 
The second cnantiomer of this chiral structure was not found. 
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TABLE III: The total energies of the most favorable clusters, per formula unit, in hartree units, at the Hartree-Fock and LDA 
level; and the same quantity for the ring structures. 



number of configu- 


HF 


LDA 


formula 


ration 






units 








1 


a 


-106.9444 


-106.6406 


2 


a 


-106.9962 


-106.6941 


3 


a 


-107.0131 


-106.7103 


1 


a 




-106.7193 




b 


-107.0187 




5 


a 




-106.7211 




d 


-107.0212 




6 


a 


-107.0279 


-106.7284 


7 


a 


-107.0281 


-106.7289 


8 


a 


-107.0315 


-106.7322 


Rings 








3 


a 


-107.0131 


-106.7103 


4 


b 


-107.0187 


-106.7152 


5 


d 


-107.0212 


-106.7172 


6 


h 


-107.0224 


-106.7184 
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FIG. 1: (Color online) Structures of (LiF)„ clusters with n=l,2,3 formula units. Lithium is displayed as a small red circle, 
fluorine as a large turquoise circle. 
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Coc, v D2/1 Dzh 

FIG. 2: (Color online) Structures of (LiF)4 clusters. For the notation, cf. fig. [T] 
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FIG. 3: (Color online) Structures of (LiF)g clusters. For the notation, cf. fig. [T] 
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FIG. 4: (Color online) Structures of (LiF)e clusters. For the notation, cf. fig. [T] 
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FIG. 5: (Color online) Structures of (LiF)7 clusters. For the notation, cf. fig. [T] 
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FIG. 6: (Color online) Structures of (LiF)s clusters. For the notation, cf. fig. [T] 
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FIG. 7: The energy (in hartree units) of the most favorable cluster for the sizes (LiF)i to (LiF)s, per formula unit (left), and 
the energy difference A/j, n = E„+i — 2E„ + E„-i (right). Note that if A/i n > 0, then two clusters of size n each are more 
favorable than two clusters with size n — 1 and n + 1. 
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FIG. 8: Tree graph for the energy barriers of the (LiF)4 cluster. 
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